Quantum rotor description of the Mott-insulator transition in the Bose-Hubbard 

model 
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We present the novel approach to the Bose-Hubbard model using the U (1) quantum rotor de- 
scription. The effective action formalism allows us to formulate a problem in the phase only action 
and obtain an analytical formulas for the critical lines. We show that the nontrivial U (1) phase field 
configurations have an impact on the phase diagrams. The topological character of the quantum 
field is governed by terms of the integer charges - winding numbers. The comparison presented 
results to recently obtained quantum Monte Carlo numerical calculations suggests that the compe- 
tition between quantum effects in strongly interacting boson systems is correctly captured by our 
model. 
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I. INTRODUCTION 

The physics of the Bose-Hubbard (BH) model was 
the subject of intensive study for a number of yearSiiSi 
More recently it has been realized that the Bose-Hubbard 
model can also be appHed to bosons trapped in so-called 
optical latticeSjiA and coarse grainingS., strong-coupling 
expansion;^ mean-field theories^ have been successfully 
applied to these systems in one-,— tviro— and three- 
dimensional lattices.— Another essentially equivalent for- 
mulation is based on the Gutzwiller wavefunction.— 
Also extension, based on a systematic strong-coupling 
approach of the Bose-Hubbard model beyond mean- 
field has been triedJ^ and experimentally confirmed.iS' 
The progress comes from better computer resources 
and more efficient algorithm allows to use the Quan- 
tum Monte Carlo (QMC) method for studies of the BH 
systems.— "i^ii^ 

An optical lattices offer remarkably clean access to a 
particular Hamiltonian and thereby serve as a model sys- 
tem for testing fundamental theoretical concepts and pro- 
viding exemplar of quantum many-body effects 1^. It is 
well known that the ground state of a system of repul- 
sively interacting bosons in a periodic potential can be ei- 
ther in a superfiuid state (SF) or in a Mott-insulting state 
(MI), characterized by integer boson densities. Because 
the phase of the order parameter and the particle num- 
ber as conjugate variables are subject to the uncertainty 
principle A(j)An ^ and the bosons can either be in 
the eigenstate of particle number or phase. The eigen- 
state of phase is a superfiuid and that of particle number 
is a localized Mott insulator. Therefore the quantum MI- 
SF phase transition takes place as the particle density is 
shifted thus facilitating emergence of the superfiuid from 
the Mott insulating state. 

The aim of this paper is to extend the mean-field 
approach for the Bose-Hubbard model in a way to in- 
clude particle number fiuctuations effects and make the 
qualitative phase diagrams in two and three dimensions 
more quantitative. Our method also improve the strong- 



coupling expansion that works well only for sufficiently 
large insulating gap.— The key point of presented ap- 
proach is to consider the representation of strongly in- 
teracting bosons as particles with attached "flux tubes". 
In a consequence a boson is the composite object. This 
introduces a conjugate U(l) phase variable, which ac- 
quires dynamic signiflcance from the boson-boson inter- 
action. To facilitate this task we employ the functional 
integral formulation of the theory that enables us to per- 
form the functional integration over fields defined on dif- 
ferent topologically equivalent classes of the U (1) group, 
i.e., with different winding numbers. An inclusion of 
the winding numbers is unavoidable in order to obtain 
a proper phase diagram. A similar method that is based 
on quantum rotor formulation was recently employed by 
one of us in the fermionic Hubbard model*^ The nice 
feature of our approach is that all the expressions and 
handling are analytic. Finally, we compare our results for 
systems at zero temperature with the outcome of the nu- 
merical simulations and found a very good agreement for 
the quantitative results regarding the behavior as we go 
from the superfiuid phase to the Mott insulating phase. 
In the framework of the introduced theory we are able to 
calculate the phase diagrams with high accuracy along 
whole critical line that separates Mott insulator - super- 
fiuid phases. Moreover, our approach gives a clear insight 
into the Bose-Hubbard model from a quantum field the- 
ory and emphasizes the impact of the topology of the 
phase variable on phase transitions. We show that the 
Coulomb interaction (as a main energy scale) is governed 
the phase transitions in the Bose-Hubbard model. 



The outline of the paper is as follows. In Sec. II we in- 
troduce the model Hamiltonian and in Sec. Ill we derive 
an effective U (1) action in the quantum rotor represen- 
tation. The aim of Sec. IV is the presentation of the 
resulting phase diagrams for two- and three-dimensional 
Bose-Hubbard systems. Finally, Section V summarizes 
our results and sets the outlook. 
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II. MODEL HAMILTONIAN 

We investigate the generic model for the Mott- 
insulator transition the Bose-Hubbard model 

H ^ ^'^''^^i'^i ^ ^) ^'^^iJ'^l^'J ^ ^^'^'^i' (1) where 

where a| and aj stands for the bosonic creation and anni- 
hilation operators that obey the canonical commutation 



relations 



ajai is the boson number 



operator on the site i, U > is the on-site repulsion 
and fi is the chemical potential which controls the num- 
ber of bosons. Here, identifies summation over the 
nearest-neighbor sites. Furthermore, Uj is the hopping 
matrix element with the dispersion for the bipartite lat- 
tice 



t (k) = 2t^coski. 



(2) 



in d dimensions. In this paper we investigate the phase 
transitions in simple cubic and square lattice. For our 
purpose we rewrite Eq. ([l]) to more suitable form 



n 



U 



(3) 



where ^i/U = ^/U +1/2 is the shifted reduced chemical 
potential. 



III. METHOD 

A. Decoupling of the Coulomb interaction 

We will adopt the method of the quantum rotor model, 
developed by one of usp2, to the BH Hamiltonian. The 
partition function of the system could be written in the 
form 



Z = y \DaDa] i 



— S\a.a\ 



(4) 



and the bosonic path-integral is taken over the complex 
fields Oi (t) with the action S given by 



Jo 



(5) 



where 



r0 Q 
Sb [a, a] = E y ^'^^ ("^^ • 



Unfortunately Hamiltonian is not quadratic in and we 
have to decouple first - the Coulomb term in Eq. ([3]) by a 



Gaussian integration over the auxiliary fields Vi (r) . The 
transformed partition function becomes 



f 


'dV 


g-52[n,V] 




2n 





(7) 



Si [a, a] 



dr 



J2 (^) 



"^Uja^ (r) aj(T) 



(8) 



and 



(9) 

After changing variables Vi (r) = Vj'^ (r) -I- the second 
part of the action takes form 



2U 



2U 



- tV^ (r) (r) 



iU 



(10) 



The field (t) could be represent as a sum of the static 
Vf (r) and periodic function V/' (t) : 



m— 1 



(11) 



where ojm = 27rm//3 (m = ±l,±2..) are the Bose- 
Matsubara frequencies. Furthermore, we introduce the 
scalar potential field which couples to the local particle 
number through the Josephson-like relation 



(12) 



where the phase field satisfies the periodicity condition 
(pi (/3) = (pi (0) as a consequence of the periodic properties 
of the V/' (r) field. We can eliminate the periodic parts 
of the fiuctuating electrochemical potential V/' (r) from 
the action replacing them by the phase field 4>i (r): 



[VaVa] e 



— Si [a, a] 



'dV^' 




. 27r 





-Sain, 4 



I, (13) 



where 



2U ^ ' ' ^tt"^ 



iU 



3 



n,(j> 



E 



dT 



(14) 



1 

2U 

i4>i (r) m (t) 



0? (r) 



iU 



4>^ (r) 



(15) 



The factor with —i dT(j)i (r) Ui (r) can be removed 
from the last equation by performing the local gauge 
transformation to the new bosonic variables as we show 
in the next subsection. 



B. Gauge transformation 

We perform the local gauge transformation to the new 
bosonic variables 







_ ai (t) 











i4>i{T) 



h (r) 



(16) 



The U(l) group governing the phase field is compact, 
i.e. 4>{t) has the topology of a circle Si, so that in- 
stanton effects can arise due to non-homotopic mappings 
of the configuration space onto the gauge group U(l). 
Therefore, we concentrate on closed paths in the imag- 
inary time (Ojl/fceT) which fall into distinct, discon- 
nected (homotopy) classes labelled by the integer winding 
numbers rii^ The chief merit of the transformation in 
Eq. (fT6| is that we have managed to cast the strongly 
correlated bosonic problem into a system of weakly in- 
teracting bosons, submerged in the bath of strongly fluc- 
tuating U(l) gauge potentials (on the high energy scale 
set by U). Now the action contains three parts: 

5i [5,5,0] = ^''dr|^5,(r)|:5,(T) 



- 5^i,,5,(r)5,(r)e-^^^^M 



(17) 



i 



2U ^ ' > iU ' 



a Jo drhdr)bdr)-- 



(18) 



i U' 



where if'ij {''') = 4>i {t) " 'Pj (''') ^-nd still we have terms 
with Vf that will be calculated in the next subsec- 
tion. Furthermore, the path-integral includes a summa- 
tion over winding numbers 



*i(0) 



(20) 



and should be performed taking phase configurations 
that satisfy boundary condition 0^ (/?) — 0^ (0) = 2imi 
where is integer. 



C. Saddle point equation 

The expectation value of the static part of the fluctu- 
ating electrochemical potential 



/[py-S]e-'S2[v^] 



(21) 



introduced in Eq. ifTTj) we calculate using the saddle point 
approximation and for U > obtain: 



i[n + U{h {r)b,{T))] 



(22) 



Now making substitution in the second part of the action 
52 [Ff] for the Vf an unique global value obtained from 
Eq. (|22l) we get flnally 



52 [5, 5] = /3^|(5,(r)6,(T 



/ drh (r) 5,; (r) . 



(23) 



The effective action is now quadratic in the bosonic vari- 
ables and can be integrated out without any difficulty 
remembering that the first term in Eq. l(23|) is simply a 
number. Therefore, the applied steepest descent method 
used to approximate integral Eq. (|2T|) allowed us to re- 
move an after effects of the auxiliary fields Vi (r) intro- 
duced in order to decouple the non-quadratic terms (in 
ai ) in Hamiltonian Eq. 



D. The partition function expressed in the phase 
fields variables 



The partition function can be expressed in form of the 
effective propagator G: 



Z = / [V<j>] 



(24) 

where exp ^— TrlnG^^^ = det G and determinant takes 
form 

det G = j [DbVb] exp I - ^ j'^ dr 



X 6,: 



d_ 



We parametrize the boson fields 

5, (r) = 5o + 5; (r) 



(25) 



(26) 
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and restrict our calculations to the phase fluctuations 
dropping the amplitude dependence. In result the inverse 
of the propagator becomes 



T - G^' (l - TGo 



(27) 



The explicit value 60 can be obtained from minimaliza- 
tion of the Hamiltonian dH (60) /dbo — where we in- 
troduced the parametrization Eq. (|26l) . Therefore, we 
write 



Go 
T 



u 



(28) 
(29) 



Kampf and Zimanyi" considered similar parametrization 
in the path-integral formulation of the coarse-graining 
procedure to the BH model. However, to obtain a criti- 
cal Hue, authors used a mean-field approach that is not 
expected to be rehable at T = and be capable to han- 
dle spatial and quantum fluctuation effects properly, es- 
pecially in two dimensions. Moreover, as we will see in 
the next sections, our results strongly depend on the di- 
mension of the system giving qualitative changing of the 
phase diagrams. 

Expanding the trace of the logarithm we have 



TrlnG^i = -Tr (in Go) - Tr (^TGo) 
1. 



-Tr 



TGo 



(30) 



with Go and T given by Eq. ([Ml and (HI). Trace over 
first term of the expansion gives us constant contribution 
to the action. From the trace over second part 

Tr (tGo) - V J^J f dT cos [4 (r) - 0, (r)] (31) 
we get the explicit form of the coefficient: 



T — hPf — 



Finally a partition function Eq. l(24|) becomes 



(32) 



(33) 



with an effective action expressed only in the phase fields 
variable 



^ph [0] 




(34) 



where = 1 if i, j are the nearest neighbors and equals 
zero otherwise. 



To proceed we replace the phase degrees of freedom 
by the complex field -01 which satisfies the quantum pe- 
riodic boundary condition V'i (/?) — i^i (0). This can be 
conveniently done using the Fadeev-Popov method with 
Dirac delta functional representation in a way used by 
Kopeo2^: 



6 (V^, - e^-^'^")) S - e-*^'(")) . (35) 



1 = 



The main idea of this approach is to attempt to gen- 
erate an effective partition function from the original 
one with cosine interaction, which incorporates the con- 
strained nature of the original variables. Thus we take tpi 
as continuous variable but constrained (on the average) 
to have the unit length: 




X e^'''^^(^ = l'^'(-)l'-^).(36) 



In Eq. I|36l) we introduced the Lagrange multipHer A 
which adds the quadratic terms (in the ipi fields) to the 
action Eq. (|34|) . Using such description is justified by 
the definition of the order parameter 



(37) 



which non-vanishing value signals a macroscopic quan- 
tum phase coherence (in our case we identify it as the SF 
state). The partition function can be written in form 



1 

where effective action ScS is given by: 

f/3 



(38) 



SeS = f drdr' [{JI,j + XS,,) S {t - t') 



Here 



7,,(t,t')]V'*V';-^A(5 (t-t'). 



(39) 



(40) 



is the two-point phase correlator associated with the or- 
der parameter field. Summarizing this part, we formu- 
lated a problem introducing an appropriate constrained 
complex order parameter field. In the next section we 
show that the presence of the nontrivial topology pos- 
sessed by the phase variable will contribute to propaga- 
tor. 
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E. Topological contribution in the correlation 
function 



The existence of the topological features of the charge 
states affects the correlation function. Because the values 
of the phases (pi which differ by 27r are equivalent thus 
we decompose the phase field in terms of a periodic field 
and term linear in r: 



(r) = if, (r) + -jriiT 



(41) 



with (P) = 0i (0) • As a result the phase correlator 
factorizes as the product of a topological term (r, r') 
depending on the integers Ui and non-topological one 



(42) 



where 



and 



(43) 



Performing the Poisson re-summation formula 



(44) 



(45) 

in (r, r') and the functional integration over the phase 
variables in (t, t') the final form of the correlator 



7y {t-,t') = 5r 



X 



-¥(».+#)%-^(".+#)(--') 



-<46) 



after Fourier transform can be written as: 

e-¥E.(n.+#)^ 



where 



(47) 



(48) 



is the partition function for the set of quantum rotors. 
The form of Eq. (|47l) assures the periodicity in the imag- 
inary time. We want to point out that another impor- 
tant property possessed by correlator Eq. (|47| comes 
out. Namely, the propagator is periodic with respect to 
/i/C/-|-l/2 which emphasizes the special role of its integer 
values. 




Figure 1: (Color online) Phase boundary between the Mott- 
insulating (MI) and superfluid (SF) phases for square (darker 
lobes) and cubic lattice in the space of the parameters t/U ~ 



IV. MOTT INSULATOR - SUPERFLUID PHASE 
TRANSITION 

The action included propagator with calculated the 
topological contribution after Fourier transform we write 

as 



(Wm) V'l 



k,m : 



k,m 



where 



rki(w,„) = A- J(k)+7-i(tj„) 



(49) 



(50) 



is the inverse of the propagator. Within the phase coher- 
ent superfiuid state the order parameter is given by 



1 



*2 = 



—y- 

k,m 



1 



J(k)+7-i(c^™)' 
where for bipartite lattices we have: 



(51) 



J(k) = 6^t(k)= (2z- + ^)t(k) 



(52) 



with the dispersion t (k) given by Eq. |[2|) and z is the lat- 
tice coordination number. The phase boundary is deter- 
mined by Eq. fSTj) from the upper limit of the eigenvalue 
spectrum max [t (k)] associated with the onset of phase 
transition. The Lagrange multiplier A "sticks" at criti- 
cality to the value Aq and stays constant in the whole 
low temperature ordered phase. The emergence of the 
critical point is signaled by the condition 



Ao - J (k = 0) + 7-1 {Ujrn=0) = 



(53) 



and by that very fact holds a converge in the constraint 
Eq. I|5ip . After summation over Matsubara frequency 
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Figure 2: (Color online) The comparison between our results 
(dashed line) and quantum Monte Carlo (black boxes). Dif- 
ference between critical values (i/t/)^;^ from our theory and 
QMC is within the range of error bars of the numerical cal- 
culations. The darker area shows the difference between the 
critical line obtained from Eq. II56P and phase boundary from 
QMC. Grey and white areas mean the Mott-insulator and 
superfluid state respectively. 



the superfluid state order parameter in the Hmit /? 
becomes 



1 - *R = 



-E 



1 



J(k=0)-./(k) 

u 



with 



"(f) 



= frac 



U 



1 

2' 



(54) 



(55) 



where frac (x) — x — [x] is the fractional part of the num- 
ber and [x] is the floor function which gives the greatest 
integer less then or equal to x. Introducing the density of 
states p (^) = A''"^ ^ K ~ ^ (k)] we obtain the critical 
line equation: 



1-^1 = 1 



+ 00 



2C(24 + # + ^)^+^M#) 



(56) 

where ^ is dimensionless parameter, ^ = Cmax— C and ^max 
stands for the maximum value of the dispersion spectrum 
t (k). The zero temperature phase diagram of the model 
calculated from Eq. ((56l) is given in Fig. [TJ We rec- 
ognize the particle-hole asymmetric - as a result of the 
model Hamiltonian Eq. ^ - Mott-insulating lobes sim- 
ilar to what was found in the literature. ^'^'^ In the MI 
phase bosons are incompressible duB/dfi ~ and local- 
ized which means that the total energy is minimized when 
each site is fllled with the same number of atoms. Increas- 
ing fluctuations in the phase system reduces fluctuations 
in the boson number on each site according to Heisen- 
berg uncertainty relation AnBA<^ > 1/2. Crossing the 
boundary Hne bosons can move from one lattice site to 
the next. The order parameter has a non- vanishing 



value and system exhibits the long-range phase coher- 
ence. This is opposite case to the Mott-insulator where 
phase coherence is lost. 

We found that our results are in great accor- 
dance with the recently published quantum Monte-Carlo 
calculations^^ (see Fig. [2]) and also improve predictions 
based on the third-order expansion mtjU that become 
inaccurate quite far from the tip?^ Furthermore, compar- 
ison with the strong coupling method^'^ indicates that it 
underestimates the critical values of tjU . For example, 
in three dimension it gives tjU — 0.029 at the tip of the 
n = 1 lobe, which is slightly lower that the value that re- 
sults form our calculation and the Monte Carlo method. 
The phase boundary is periodic with respect io n/U with 
fixed integer fiUing depending on the value of the chem- 
ical potential p,. The vicinity of the lobe tip (t/J/)^.^.;^, 
corresponding to the MI-SF transition in the commen- 
surate system is shifting from value 0.4 to 0.5 when we 
change a dimension of the lattice from three- (3Z?) to two- 
dimensional (2D). Moreover, we see that the qualitative 
shape of the lobes is not the same for 2D and 3Z) cases 
and steeper for the two-dimensional system. Analysis of 
the one-dimensional systems is not possible in presented 
approach because for dimensions d < 2 it does not ex- 
hibit the phase transition at finite temperatures T > 0, 
in agreement with the Mermin- Wagner theorem.'^ 

Finally a comment regarding the critical behavior of 
the model in our quantum rotor approach is in order. 
To extract the near-critical form of the propagator it is 
this sufficient to perform an expansion in terms of the 
momentum k and frequency oJm in Eq. l(50|) . In the 
T ^ limit, with the help of Eq. l(47|) . after proper 
re-scaling of the fields ■(/'k.tJm we find 



r; 



; ^ (t^m) = r + k^ + LU„i + iujjn + V (^jj^ 



(57) 



Here, r ^ J (k = 0) — A is the critical "mass" parameter 
that vanishes at the phase transition boundary, k^ = /c • /c 
while V {%) is given by Eq. (|55l) . Due to the quantum 
nature of the problem, the scaling of the spatial degrees 
of freedom k ^ k' = sk implies the scaling for frequen- 
cies in a form uJm — + = s^ujm with the dynamical 
critical exponent z. At the tips of the lobes in the t/U- 
p,/U phase diagram in Fig. [2] one has v{fi/U) = 0, so 
that T^^{uJm) ~ fc^ +i^m with space-time isotropy giving 
z — 1. However, the other points on the critical line with 
non-vanishing v [fJ-lU) refiect the absence of the particle- 
hole symmetry due to the imaginary term involving iojm- 
In this case the higher order term WmlV'k.Wm P becomes 
irrelevant and can be ignored, while the critical form of 
the propagator (|57|) reads r^^(w„i) ^ + iv {n/U)uJm- 
Now, the scaHng requires z = 2 as a result of the 
momentum-frequency anisotropy. 



7 




Figure 3: (Color online) Boson occupation number ub at 
T = for three-dimensional simple cubic lattice in the space 
of parameters - chemical potential fi/U and hopping t/U. The 
Mott insulator is found within each lobe of integer boson den- 
sity. Inside the first lob on the left the occupation number ub 
is equal one, two and three in second and third step respec- 
tively. 



A. Boson occupation number 

The effects of the fixed boson number in the system 
defined by 



(58) 



are included in our theory because of the source term con- 
taining chemical potential l^'Ylii dTbi{T)hi{T) in ac- 
tion Eq. (|23l) . By differentiating the partition function 
Eq. lfT3|) (after carrying out a gauge transformation and 
change the variable in the action) we obtain 



(fli (r) ai (r)) = 



d\uZ 
dp. 



1 



{vf) + U.)-ifi 



(59) 

Inserting in above a static part of the electrochemical 
potential Eq. l(22|) we find the boson density 



riB 



N ^ 



{k (r) h (r)) + ^ (cf>, 



(60) 



When the phase stiffness vanishes J — the bosonic 
contribution to the free energy is given by 



1 

']3N 



In / [V(f>,] e 



which is simply the contribution from the "free" rotor ac- 
tion. Now, we again decompose the phase field in terms 
of a periodic field and term linear in r. Calculating in- 
tegral Eq. ((6T1) we get in the limit T ^ an analytical 
solution 



(m)|j=o = 



U 



from which we recognize a steps of fixed integer filling of 
bosons (set t/U = in the Fig. [S]). 

The calculations of a phase diagram for interaction 
problem t/U ^ are more complicated since spatial cor- 
relations have to be included, as well. However our model 
is expressed in terms of the complex field Eq. ([351 which 
is now very helpful. The result for the boson density tib 
within the region of superfuidity is given by the expres- 
sion 



rtB 



u 



1 



- 2^- 



B' 



(63) 



(61) 



where non- vanishing value of the order parameter is 
calculated from Eq. ((56l) . We see in the Fig. [3] that the 
competition between kinetic and interaction energy is the 
foundations of the quantum phase transitions in the BH 
model. Increasing the value of the hopping term (reduc- 
ing the interaction energy) leads to delocalization of the 
bosons thus the sharp steps of the MI state become in- 
distinct and in consequence system is superfiuid. In Fig. 
[3] we observe the appearance of the Mott-insulating lobes 
corresponding to curves from Fig. [TJ The MI has a gap 
to density excitations and is an incompressible (density 
plateaus in Fig. ^ thus the chemical potential can be 
changed within a gap without changing the density. At 
the tip of the lobe at fixed integer density the transi- 
tion is driven by the change of the t/U ratio in a system 
composed of a fixed number of bosons. Such a transi- 
tion in a d-dimensional BH model lies in the universality 
class of the {d + l)-dimensional XY spin model. Remain 
possibilities that the system can cross the superfiuid - 
Mott insulator phase boundary are called generic (when 
we add / subtract a small number of particles) and do not 
belong with the universality class of the XY spin model, 
so are characterized by different critical exponents. 

The possibility to describe both the Mott and SF 
phases in two dimensions can be also done using the 
strong-coupling expansion method.— The obtained re- 
sults are qualitatively comparable to our phase diagrams 
and show that both extensions of the Bose-Hubbard 
model beyond mean-field succeed in catching strongly 
interacting systems. Besides, authors calculated the ex- 
citation energies and spectral weight and provided ana- 
lytical formulas which expanded can be useful to deter- 
mine the expected second order term in the momentum 
distribution. •!£ Furthermore it seems that presented ap- 
proaches can be in principle applied to more complicated 
situations. 
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V. SUMMARY AND OUTLOOK 

In this paper we have presented a study of the Mott- 
insulator transition of the Bose-Hubbard model. To ana- 
lyze a quantum phase transitions beyond mean-field the- 
ory we employed a U(l) quantum rotor approach and 
a path-integral formulation of quantum mechanics in- 
cluding a summation over a topological charge, explic- 
itly tailored for the BH Hamiltonian. The effective ac- 
tion formalism allows us to formulate a problem in the 
phase only action and obtain an analytical formulas for 
the critical lines. We have compared obtained results 
to existing numerical calculations and found them in a 
very good agreement. The formalism adopted here can 
be extended and apphed to the other systems systemat- 



ically. Especially the effect of the competition between 
quantum effects in finite temperatures focuses our atten- 
tion. Considerations different geometries of the lattices 
are possible in the frame of our approach, as well. These 
topics will be considered in future publications. 
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